Propagation of dark soliton interacting with domain wall in two immiscible Bose–Einstein condensates*

Project supported by the National Natural Science Foundation of China (Grant Nos. 61565007, 11875149, 11747079, and 11874127), the Science Fund from the Department of Science and Technology of Jiangxi Province, China (Grant Nos. 20162BCB23049 and 20171ACB21045), the Youth Jinggang Scholars Program in Jiangxi Province, China, and the Program of Qingjiang Excellent Yong Talents, Jiangxi University of Science and Technology, China.

Zheng Lang1, 2, Zhang Yi-Cai3, Liu Chao-Fei2, †
Jiangxi Provincial Education Examination Authority, Nanchang 330038, China
School of Science, Jiangxi University of Science and Technology, Ganzhou 341000, China
School of Physics and Electronic Engineering, Guangzhou University, Guangzhou 510006, China

 

† Corresponding author. E-mail: liuchaofei0809@163.com

Project supported by the National Natural Science Foundation of China (Grant Nos. 61565007, 11875149, 11747079, and 11874127), the Science Fund from the Department of Science and Technology of Jiangxi Province, China (Grant Nos. 20162BCB23049 and 20171ACB21045), the Youth Jinggang Scholars Program in Jiangxi Province, China, and the Program of Qingjiang Excellent Yong Talents, Jiangxi University of Science and Technology, China.

Abstract

Reflection and transmission are two behaviors of wave propagating to an interface. The immiscible binary mixtures of Bose–Einstein condensates can form the symmetry-breaking state, in which the domain wall on the center can serve as the interface. In this study, we explore in detail the propagation of a dark soliton interacting with the domain wall in the harmonic trap. We find that the low-energy dark soliton is easy to form the transmission and the high-energy dark soliton trends to reflect from the domain wall. Both reflection and transmission of dark soliton on the domain wall induce the sound radiation. But the sound radiation in the reflection derives from the collective oscillation of condensates, and it in the transmission comes not only from the collective oscillation, but also from the condensate filling in the dark soliton.

1. Introduction

Transmission refers to the part of the wave that can continue to propagate through the interface when the wave propagates to the interface between two media. Usually, when the wave encounters an interface, part of it is reflected back from the interface, and the remaining part goes through it into another medium. What happens when the system is a superfluid medium and a dark soliton propagates, and the interface is a domain wall composed of two incompatible superfluid media? Dark soliton, as a fundamental excitation in the repulsive atomic Bose–Einstein condensate (BEC), has been studied extensively and attracted more and more attention.[18] Dark soliton is a localized one-dimensional (1D) defect characterized by a notch in the ambient condensate density. As is well known, dark soliton, accompanied by a phase jump, results from a balance between the defocusing dispersion and the focusing repulsive nonlinear interaction. Thus, dark soliton can preserve its localized shape during propagation and display particle-like characters.[912]

Almost all investigations about dark soliton concentrate on the single-component BEC. Experimentally, the dark soliton has been generated in diluted BEC by the phase-imprinting method[1,2] and the perturbation of density.[5] In the phase-imprinting method, part of the BEC is irradiated by a far detuned laser for a short time, which results in the phase shift without sensitive density fluctuation. In this experiment, not only one dark soliton but also some sound waves are created. After that, due to the inherent instability and the transversal excitation, the dark soliton has been observed to degenerate into vortex rings. In theory, the BEC formed by diluted atoms can be described by Gross–Pitaevskii (GP) equation, which is a typical nonlinear system. When the external potential is zero, the GP equation supports the solution of the dark soliton in one dimension. This is the reason that a large number of previous studies of dark solitons prefer to consider the system of GP equation. Meanwhile, the studies about dark solitons mainly focus on some new solutions and the exploration of the instability. In 2001, Öhberg and Santos extended the investigation of dark soliton on single-component BECs to the two-component case.[13] They have considered two scenarios: the miscible components and the immiscible case. For miscible components, they analyzed the creation of solitons in different components via phase-imprinting, and the soliton–soliton interaction in different situations, which include the formation of a soliton–soliton binary system and free single solitons accompanied by a filling of the other component. For immiscible components, they showed that a dark soliton can be transferred from one component to another at the domain wall when exceeding a critical velocity.

In fact, binary mixtures of BECs possess complex spatial structures.[1421] To determine the density profile, Ho and Shenoy have first presented a simple algorithm within the Thomas–Fermi (TF) approximation.[14] More accurately, Hartree–Fock theory is provided to achieve the rich spatial structures.[1517] In general, the interspecies interaction coefficient g12 plays an important role in determining the structure of ground state. When the inequality is satisfied, the two BECs are miscible; and when , the trapped BECs are immiscible due to the strong interspecies repulsion. Furthermore, the spatial patterns of the immiscible mixtures have two basic configurations: the symmetry-breaking state (SBS) and the symmetry-preserving state (SPS). And one of the configurations can be viewed as the metastable state.[2224] In Ref. [13], Öhberg and Santos have used a box trap where the BECs form the SPS, which display the dark soliton propagation from one component to another.

In this study, we choose the SBS configuration as a platform to demonstrate the propagation of dark soliton. The separated BECs are confined in a harmonic trap and the dark soliton is initially created off-centre. In this situation, the domain wall will be located near the center of the trap. When the soliton is far from the domain wall, it will act as a single-component BEC.[25,26] Hence, the soliton accelerates towards the trap center, and it gains a big velocity to pass through the domain wall. These are the main differences from the model in Ref. [13]. Furthermore, we will examine in detail the evolution of the soliton energy and the soliton depth. Especially, we find the soliton carries some atoms of the other component into the BEC as they pass through the domain wall. This leads the the soliton to dissipate eventually.

2. Model and basic equations

One would consider the binary mixtures in a three-dimensional (3D) trap potential which can be described by

where mi is the atomic mass, and ωi and ωi are the longitudinal and transverse trapping frequencies, respectively. If the longitudinal trapping frequency is much smaller than the transverse trapping frequency, the trapping potential is cigar-shaped. In addition, ωi does not affect the transverse component of the wave functions when the two-body interaction energy is much smaller than ℏω. These features allow us to investigate the propagation of soliton in the 1D space.

For sufficiently low temperatures, the system can be described by the two coupled Gross–Pitaevskii (GP) equations[1321]

Here, ψ1 and ψ2 denote the macroscopic wave functions of the two BECs respectively. The intraspecies interaction coefficients are g1 = 4πℏ2a1/m1 and g2 = 4πℏ2a2/m2, where a1 and a2 are the scattering length of the species 1 and the species 2 respectively; and the interspecies interaction coefficient is g12 = 2πℏ2a12/[m1m2/(m1 + m2)], where a12 is the scattering length between species 1 and species 2. In addition, each wave function is normalized by the number of particles Ni as ∫ dx|ψi(x)|2 = Ni.

Because the two components construct the SBS, we can set a dark solution only to be in one of the BECs (see the scheme of our model in Fig. 1). Usually, dark soliton solutions are supported by the 1D single-component GP equation in the absence of external confinement.[12,27] Here, we assume that the dark soliton is in BEC 1. Therefore, the dark soliton with speed v, position (xvt) and a uniform background density n1 has the form,

where , is the chemical potential, is the Bogoliubov speed of sound, and is the healing length characterizing the size of the soliton.

Fig. 1. Initial condensate densities in harmonic trap (not shown) with a stationary soliton being located at x0 = −30ξ, length unit being . Corresponding unit in the following pictures is the same as in this picture. When dark solitons are not considered, density distribution of two BECs (short as BEC 1 and BEC 2) comes from imaginary-time propagation.

In numerical experiments, the soliton energy can be calculated by integrating the two coupled GP energy functions

across the soliton region (X ± 5ξ, where X is the instantaneous position of the local density minimum) and subtracting the corresponding contribution of the background fluids.

3. Initial condition

Because there are two different components in our model, the initial profile cannot be easily obtained in the TF approximation. Thus, we have to firstly obtain a stable SBS configuration. It can be achieved by mixing two initially separated BECs into one trap.[24] In Ref. [24], we have confirmed that the direct mixture process will firstly form the SBS. Therefore, we use the mixture model in Ref. [24] with the external potential , i = 1,2. Unloading the splitting potential δ(x), and propagating the corresponding TF wave function in imaginary time, we can obtain the smoothed density profile n1(x) and n2(x) which are relatively stable in the harmonic trap (see Fig. 1).

For simplicity, we assume that the confining potentials for the two components are equal, i.e., . Here, we set for the external potentials. Meanwhile, the interaction coefficients are g1 = 1, g2 = 0.8g1, and g12 = 1.2g1. We use = 1 and the atomic mass m1 = m2 = m = 1. In this situation, the spatial extent of the system is characterized by the healing length , and the time unit is ξ/c. Note that we choose μ1 = 1 and obtain μ2 = 0.8618μ1 by calculation to satisfy N1 = N2.

We now estimate the parameters for a realistic experiment. We assume that BEC 1 is formed by 23Na with m = 38.18 × 10−27 kg and a = 2.8 nm.[28,29] The tight transverse confining frequency is ω = 5 000 × 2π Hz and the 1D peak condensate density is n1D = 108 m−1. Thus, the longitudinal confining frequency is 39.7 × 2π Hz. Our space and time units correspond to 0.4 μm and 5.7 × 10−5 s respectively. The two BECs have the same number of atoms N1N2 ≈ 2.6 × 104.

4. Results
4.1. Reflection oscillations of dark solitons interacting with domain wall

Usually, the dark soliton deviating from the position of the lowest potential energy will oscillate in the harmonic potential well. It behaviors like the oscillating motion of a small ball away from the equilibrium position in the harmonic trap. Here, there is a domain wall formed by the two phase-separated BECs near the the position of the lowest potential energy the trap, and the dark soliton will inevitably encounter and interact with the domain wall in its motion.

When the initial size of the dark soliton is not particularly large to deviate from the domain wall, the speed of the dark soliton is very small when it reaches the domain wall. This inevitably leads to the fact that although it can interact with the domain wall, it cannot reach the BEC 2 by passing the domain wall. Figure 2 shows this situation, where we initially set the dark solitons to be at the position of −6.1ξ. Figure 2(a) is the density of the BEC 1, in which dark solitons oscillate slightly. Meanwhile, in the process of dark soliton oscillation, BEC 1 keeps its original shape basically unchanged. Figure 2(b) shows the density of the BEC 2, which is not affected by the dark solitons. Figure 2(c) displays the total density of the BECs. Figure 2(d) presents the result from which we have removed the background of the BEC. Here, we can clearly find the sound radiation in both the BEC 1 and BEC 2.

Fig. 2. Reflection oscillations of dark solitons interacting with domain wall, where initial position of dark soliton is set to be at −6.1ξ. (a) Density of condensate 1, (b) density of condensate 2, (c) density of whole condensates, and (d) sound radiation of condensate.

It is not surprising to obtain the sound radiation in BEC 1 because the dark soliton would oscillate in it. However, dark solitons do not propagate in the BEC 2. Why can the sound radiation appear in the BEC 2?

Figure 3(a) shows the oscillation trajectories of the dark solitons. It can be seen that the dark soliton only oscillates below the x axis. In fact, this kind of oscillation is different from the normal one where the dark soliton moves back and forth from positive coordinate to negative coordinate. At the same time, the oscillation is also similar to the sinusoidal trajectory. Figure 3(b) shows the variation of the energy of the dark soliton. The interaction between dark soliton and sound wave makes its energy change slightly and periodically. Figure 3(c) shows the oscillation of the two phase-separated BECs. The blue curve shows the oscillation of BEC 1, which has a translation of 73.6ξ. Since the dark soliton moves in BEC 1, it is easy to understand the centroid oscillation of BEC 1, which is always opposite to the dark soliton’s oscillation. We can also see that the centroid oscillation of BEC 2. Meanwhile, the centroid oscillation of the two BECs are almost synchronous.

Fig. 3. Properties of reflection oscillation of dark solitons interacting with domain wall. (a) Oscillation route of dark soliton, (b) energy of dark soliton, and (c) oscillation of two-component condensate.

Now, we can explain why there is sound radiation in BEC 2. The dark soliton oscillates in the BEC 1, leading to the centroid oscillation of the BEC 1. At the same time, the centroid oscillation of the BEC 1 also induces the corresponding centroid oscillation of BEC 2. During this oscillation, the density modulation of BEC 2 makes the sound radiation occur in the harmonic trap.

4.2. Transmission propagation of dark solitons through domain wall

When the initial position of the dark soliton is far from the domain wall, it can move very fast when it reaches the domain wall. At this time, we can see the phenomenon that the dark soliton passes through the domain wall and interacts with it. Here, we set the dark soliton of the initial velocity to be 0 at −30. Figure 4 shows the long-time oscillation behaviors. Figure 4(a) is the density of BEC 1, in which dark solitons oscillate with big amplitude. Obviously, dark soliton transmits into BEC 2 [see Fig. 4(b)] and oscillates back and forth in both BECs. After several cycles of oscillation, we can also see some white trajectories of dark soliton in BEC 1 within the region of X > 0. At X < 0, the trajectories of dark soliton also appear in BEC 2. In Fig. 4(c), we can see the total trajectory of dark solitons and the shift of the domain wall near the trap center. Since dark soliton oscillates in both BECs, in Fig. 4(d) the sound radiation is obvious in both BECs after removing the background of the BECs.

Fig. 4. Transmission oscillation of dark soliton interacting with domain wall, where initial position of the dark soliton is set to be at −30ξ and the initial velocity of the dark soliton is set to be 0. (a) Density of BEC 1, (b) density of BEC 2, (c) density of the whole condensates, and (d) sound radiation of condensates.

Figure 5(a) shows the trajectory of dark soliton in the two BECs. The amplitude of oscillation is not particularly fixed. Figure 5(b) shows the variation of the dark soliton energy. There are many platforms with varying heights and some sagging lines. In fact, the platform corresponds to the energy of the dark soliton when it is completely in one BEC, and the sagging line is the result when the dark soliton passes through the domain wall. Figure 5(c) shows the variation of dark soliton depth. There are also many platforms and some sagging lines in these results. The platform corresponds to the depth when the dark soliton is completely in a certain BEC, and the sagging lines are the result of the depth of the dark soliton when it passes through the domain wall. Figure 5(d) is the centroid oscillation of the two BECs. The red curve is the result of BEC 1 (Note that we obtain the results by a translation of 73.6ξ for comparing with those of BEC 2). Both BEC1 and BEC 2 oscillate as the dark soliton moves in them. Meanwhile, the collective oscillations of the two BECs are basically synchronous.

Fig. 5. Properties of transmission oscillations of dark solitons interacting with domain wall. (a) Oscillation trajectory of dark soliton, (b) energy of dark soliton, (c) depth of dark soliton, and (d) oscillation of two-component condensate.

Figure 6 shows the enlarged plots about the propagation of the dark soliton. The dark soliton is initially created off-center in the trap. From its oscillating trajectory, we can see that the dark soliton in this system can display particle-like characters, just like those in the single-component BEC. Furthermore, the domain wall is not fixed. When the soliton first passes through the domain wall, this position changes suddenly, bounding from 1.7ξ to −0.4ξ. It is just the so-called back-action in Ref. [13]. Certainly, the BECs oscillate due to the motion of the soliton. This causes the domain wall to shift slightly.

Fig. 6. Propagation of dark soliton in two immiscible BECs. The initial position of the soliton is at x0 = −30ξ and the initial speed is zero. Time unit is ξ/c. Corresponding unit in the following figures is the same as in this picture.

To fully display the evolution of dark soliton, we also examine the soliton velocity [Fig. 7(a)], the soliton depth [Fig. 7(b)] and the soliton energy [Fig. 7(c)]. Initially, the soliton experiences a force: the force accelerates the soliton towards the trap centre. Its velocity becomes complex as the soliton encounters the domain wall. After it passes through the domain wall, the soliton enters into the BEC 2 and becomes decelerated until its moving direction changes [see Fig. 7(a)]. And then, soliton returns from BEC 2. Therefore, the dark soliton is expected to oscillate back and forth in the trap.

Fig. 7. (a) Velocity of soliton versus time, (b) depth of soliton versus time, and (c) soliton energy versus time. Initial position of the soliton is at x0 = −30ξ and initial speed of the soliton is zero. Energy is in units of μ1. Corresponding unit in the following figures is the same as in this figure.

In fact, the depth of the soliton density depression stays constant when the soliton is singly embedded in BEC 1 or BEC 2 (see Fig. 7(b)). Note that the soliton depth is the difference between the corresponding background BECs and the local density is a minimum value. Especially, the intraspecies interaction strength of BEC 2 is less than that of BEC 1. Hence, the peak density of the BEC 2 is bigger than that of the BEC 1. Accordingly, the soliton in BEC 2 produces a bigger depth [see Fig. 7(b)]. Similarly, the soliton energy is roughly conserved when it is completely in BEC 1 or BEC 2. The difference between values about soliton energy in the two components results from the different intraspecies interaction strength. Clearly, the intraspecies interaction strength of BEC 2 is smaller than that of BEC 1. Therefore, the soliton energy is smaller to that of in BEC 2.

Generally speaking, the change of soliton is unusual when it encounters the domain wall. Note the trapped BECs are inhomogeneous, so the initial soliton becomes asymmetrically deformed and tries to adjust to the inhomogeneous background by radiating counter propagating sound pulses.[30] Thus, the soliton energy decreases weakly in BEC 1 at first. Here, the value seems to be conserved [see Fig. 7(c)]. When soliton approaches to the domain wall (t = 150ξ/c), the soliton energy increases because the addition of the strong interspecies interaction. When the soliton is located in the domain wall region, the soliton depth suddenly decreases to a minimum value due to the strong interspecies interaction. Accordingly, the soliton energy decreases to a minimum value. Until the soliton transmits into the BEC 2, its energy attempts to revert to a new level and tends to be conserved. The behaviors of soliton transferring from BEC 2 into BEC 1 are similar. Furthermore, the soliton depth has a small decrease instead of regaining the initial value as the soliton returns to BEC 1 [see Fig. 7(b)].

Now, we further focus on the detailed process that the soliton propagates from BEC 1 into BEC 2 (Fig. 8). In Ref. [13], Öhberg and Santos have shown that a dark soliton can be transferred from one component to another at the domain wall when it exceeds a critical velocity. In our model, the domain wall is located near the origin. The soliton velocity is almost the maximum value at this position. So it easily passes through the domain wall. A very interesting point is that the dark soliton can carry some atoms of BEC 1 into BEC 2 [see Figs. 8(e) and 8(f)]. Similarly, as the soliton propagates from BEC 2 into BEC 1, it carries some atoms of BEC 2 into BEC 1.

Fig. 8. Soliton propagation from BEC 1 (black line) into BEC 2 (red line) at t = 100ξ/c (a), 150ξ/c (b), 180ξ/c (c), 250ξ/c (d); (e) enlarged figure of panels (c), (f) enlarged figure of panel (d).

For the soliton in single-component BEC trapped in harmonic potential, the oscillation behavior is robust if the system is ideal. If there is perturbation such as optical lattices,[31] damping,[32] and finite temperature,[33,34] the oscillation amplitude tends to increase regularly as the soliton energy decreases. Meanwhile, the soliton depth decreases until the soliton disappears. Therefore, the oscillation behavior is strange (see Figs. 4 and 5). The oscillation amplitude does not increase regularly. Reversely, we can see that the amplitude of soliton oscillation is not stable and the oscillation period tends to increase (see Fig. 5). Certainly, when the soliton becomes weak enough, it is unable to pass through the domain wall. The dark soliton disappears eventually. Therefore, the SBS configuration is not an ideal system for soliton oscillation, although we do not add any dissipation.

In the above text, we have revealed in detail the propagation of soliton from one BEC into the other. Now, we can explain the soliton behaviors in the long time evolution. In Figs. 8(e) and 8(f), we note that the soliton carries some atoms of the other components into the embedded BEC. This is a key point for causing the above phenomena to happen. Figure 5(b) indicates that the soliton energy tends to decrease in the long time evolution. Figure 5(c) shows that the corresponding soliton depth tends to decrease. The soliton carries more and more atoms of BEC 1 (2) into BEC 2 (1). Undoubtedly, the added atoms will perturb the soliton and the soliton cannot function normally. Because the soliton carries more and more atoms of the other components, the soliton depth tends to decrease. Meanwhile, we can imagine that the added atoms would cause the soliton to move difficultly. Hence, the oscillation amplitude of soliton does not approach to the edge of the trap regularly as the soliton energy tends to decrease.

We have taken the interaction coefficients g1 = 1, g2 = 0.8g1, and g12 = 1.2g1 for example. We study the propagation of soliton in two immiscible BECs, which form the SBS configuration. Unlike the scenario in Ref. [13], our investigation concentrates on the change of soliton energy and soliton depth. Especially for the soliton energy, this work first indicates its evolution in the binary mixture. In fact, the propagation properties of dark soliton are also dependent on other factors such as the trapping potential, the atom number of the two components, the interspecies interaction and the intraspecies interaction. If the interspecies interaction is large enough, the soliton would be very hard to transfer from one component to another. If the two BECs have different numbers of atoms, the two BECs will form an extremely asymmetric structure. Even the domain wall is far from the trap center, it will be difficult for the soliton to reach and pass through.

4.3. Critical transmission of dark solitons interacting with domain wall

In the above study, we gave two kinds of interactions between dark soliton and the domain wall. One is the reflection and the other is the transmission. It is also very easy to realize a critical case of transmission for the dark soliton interacting with the domain wall. Now, there come some problems. When will it occur? Where should the initial position of the static dark soliton be set to be? What happens in this situation?

In our simulation, we find that the transmission of dark solitons cannot occur at 6.1ξ. When the initial position of dark solitons is 6.2ξ, we find the critical transmission will occur. Figure 9 shows the result of critical transmission. Like Figs. 2 and 4, figures 9(a) and 9(b) show the density evolution of BEC 1 and BEC 2, respectively. The density of the whole BECs is shown in Fig. 9(c). And figure 9(d) shows the sound radiation from which we have removed the background BECs. In fact, we can easily find that the trajectories of dark solitons in the two BECs are different. In Fig. 9(a), the oscillation amplitude of dark solitons is smaller than that in Fig. 9(b).

Fig. 9. Critical transmission of dark solitons interacting with domain wall, where initial position of dark soliton is set to be at −6.2ξ. (a) Density of BEC 1, (b) density of BEC 2, (c) density of the whole condensates, and (d) sound radiation of condensates.

Figure 10(a) shows the trajectory of the dark soliton. The oscillation amplitude of dark soliton in BEC 1 is about 6.2ξ but in BEC 2 the oscillation amplitude is close to 20ξ. Figure 10(a) clearly shows the difference in oscillation amplitude. The energy of dark solitons in different BECs will change abruptly. Obviously, when oscillation amplitude of the dark soliton is larger, its energy is lower (see Fig. 10(b)). On the contrary, the oscillation amplitude is smaller and the energy is larger. Even if it moves into the two different BECs, the inverse relationship between energy and oscillation amplitude is obvious. Figure 10(c) shows the changes of the centroid position of the two BECs.

Fig. 10. Properties of critical transmission of dark soliton interacting with domain wall: (a) oscillation route of dark soliton, (b) energy of dark soliton, and (c) oscillation of the two-component condensates.

Snake instability is a dynamical behavior of dark solitons to decay into vortices in quasi-one-dimensional (quasi-1D), 3D and some real systems.[3537] Firstly, the solution of dark soliton in GP equation for BEC usually suits the homogeneous system,[30,32] that is, there is no external potential and the density of the BEC is uniform everywhere. Only in this situation there is the exact mathematical solution of dark soliton. But, the real system is very difficult to achieve this condition. In addition, the exact solution of the dark soliton would be used as the initial condition when simulating the dark soliton in the quasi-1D, 3D, and some real systems. At this time, the dark soliton behaves no longer in the ideal 1D homogeneous BEC.[38,39] Now we take the quasi-1D BEC for example. The dark soliton has initially a straight snake-shape perpendicular to the flat (quasi-1D) condensate. The quasi-1D dark soliton is unstable, so the straight snake-shape is twisted snake-shaped. Finally, the twisted snake-shaped dark soliton would decay into vortices. This is how the snake instability occurs in the evolution of dark soliton. In this work, we consider the 1D BEC rather than the quasi-1D BEC, so the snake instability does not occur in the present study.[38,39]

5. Conclusions

In this study, we investigated the propagation of dark soliton in two immiscible BECs, which form the SBS configuration in the harmonic trap. The soliton is initially off-center and accelerated to the trap center. When the static dark soliton is initially located below some critical position, it reflects from the domain wall and periodically oscillates only in one component of the BEC. But the sound rediation would occur in both condensates due to the collective oscillation. When the dark soliton is located above the critical position, it passes through the domain wall to oscillate between the two condensates. We find that the dark soliton can carry some atoms of the other component. This gives rise to the sound radiation and the decrease of the solion depth and the loss of the soliton energy.

Reference
[1] Burger S Bongs K Dettmer S Ertmer W Sengstock K Sanpera A Shlyapnikov G V Lewenstein M 1999 Phys. Rev. Lett. 83 5198
[2] Denschlag J Simsarian J E Feder D L Clark C W Collins L A Cubizolles J Deng L Hagley E W Helmerson K Reinhardt W P Rolston S L Schneider B I Phillips W D 2000 Science 287 97
[3] Anderson B P Haljan P C Regal C A Feder D L Collins L A Clark C W Cornell E A 2001 Phys. Rev. Lett. 86 2926
[4] Engels P Atherton C 2007 Phys. Rev. Lett. 99 160405
[5] Dutton Z Budde M Slowe C Hau L V 2001 Science 293 663
[6] Carr L D Br J Burger S Sanpera A 2001 Phys. Rev. A 63 051601(R)
[7] Brazhnyi V A Kamchatnov A M 2003 Phys. Rev. A 68 043614
[8] Burger S Carr L D Öhberg P Sengstock K Sanpera A 2002 Phys. Rev. A 65 043611
[9] Reinhardt W P Clark C W 1997 J. Phys. B: At. Mol. Opt. Phys. 30 L785
[10] Busch T Anglin J R 2000 Phys. Rev. Lett. 84 2298
[11] Pelinovsky D E Frantzeskakis D J Kevrekidis P G 2005 Phys. Rev. E 72 016615
[12] Brazhnyi V A Konotop V V Pitaevskii L P 2006 Phys. Rev. A 73 053601
[13] Öhberg P Santos L 2001 Phys. Rev. Lett. 86 2918
[14] Ho T L Shenoy V B 1996 Phys. Rev. Lett. 77 3276
[15] Esry B D Greene C H Burke J P Jr Bohn J L 1997 Phys. Rev. Lett. 78 3594
[16] Timmermans E 1998 Phys. Rev. Lett. 81 5718
[17] Öhberg P Stenholm S 1998 Phys. Rev. A 57 1272
[18] Pu H Bigelow N P 1998 Phys. Rev. Lett. 80 1130
[19] Graham R Walls D 1998 Phys. Rev. A 57 484
[20] Chui S T Ao P 1999 Phys. Rev. A 59 1473
[21] Ao P Chui S T 1998 Phys. Rev. A 58 4836
[22] Öhberg P 1999 Phys. Rev. A 59 634
[23] Kasamatsu K Yasui Y Tsubota M 2001 Phys. Rev. A 64 053605
[24] Liu C F Tang Y 2009 Eur. Phys. J. B 70 193
[25] Liu C F Hu K Hu T Tang Y 2011 Chin. Phys. B 20 010309
[26] Liu C F Pan X Q Zhang G Y 2016 Journal of Jiangxi University of Science and Technology 37 102 in Chinese
[27] Konotop V V Pitaevskii L 2004 Phys. Rev. Lett. 93 240403
[28] Liu C F Liu W M 2012 Phys. Rev. A 86 033602
[29] Samuelis C Tiesinga E Laue T Elbs M Knöckel H Tiemann E 2000 Phys. Rev. A 63 012710
[30] Parker N G Proukakis N P Leadbeater M Adams C S 2003 Phys. Rev. Lett. 90 220401
[31] Parker N G Proukakis N P Barenghi C F Adams C S 2004 J. Phys. B: At. Mol. Opt. Phys. 37 S175
[32] Proukakis N P Parker N G Barenghi C F Adams C S 2004 Phys. Rev. Lett. 93 130408
[33] Jackson B Proukakis N P Barenghi C F 2007 Phys. Rev. A 75 051601(R)
[34] Liu C F Fan H Zhang Y C Wang D S Liu W M 2012 Phys. Rev. A 86 053616
[35] Feder D L Pindzola M S Collins L A Schneider B I Clark C W 2000 Phys. Rev. A 62 053606
[36] Huang G Makarov V A Velarde M G 2003 Phys. Rev. A 67 023604
[37] Br J Reinhardt W P 2002 Phys. Rev. A 65 043612
[38] Liu C F Lu M Liu W Q 2012 Phys. Lett. A 376 188
[39] Liu C F Hu K Hu T Tang Y 2010 Phys. Lett. A 374 2089